Optimal transport on supply-demand networks 
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Previously, transport networks are usually treated as homogeneous networks, that is, every node 
has the same function, simultaneously providing and requiring resources. However, some real net- 
works, such as power grid and supply chain networks, show a far different scenario in which the 
nodes are classified into two categories: the supply nodes provide some kinds of services, while 
the demand nodes require them. In this paper, we propose a general transport model for those 
supply-demand networks, associated with a criterion to quantify their transport capacities. In a 
supply-demand network with heterogenous degree distribution, its transport capacity strongly de- 
pends on the locations of supply nodes. We therefore design a simulated annealing algorithm to find 
the optimal configuration of supply nodes, which remarkably enhances the transport capacity, and 
outperforms the degree target algorithm, the betweenness target algorithm, and the greedy method. 
This work provides a start point for systematically analyzing and optimizing transport dynamics 
on supply-demand networks. 

PACS numbers: 89.75.Hc, 05.60.-k, 89.20.Hh 



I. INTRODUCTION 

Network transport has attracted increasing attention 
in recent years (see the review articles [ll, Q and the refer- 
ences therein) . Indeed, it describes a large number of nat- 
ural phenomena and technological processes, such as sub- 
stance flow in a metabolic network, power transmission in 
an electric network, information propagation in the Inter- 
net, and so on. A matter of prime importance is to make 
the transport processes more effective and efficient, cor- 
responding to maximizing the global capacity and min- 
imizing the average delivery time. Previous works ad- 
dressed this issue can be roughly classified into two cat- 
egories: one concerns the optimization/modification of 
underlying topology [1, 0, [Bl] , while the other focuses on 
the design of highly efScient transport/routing protocols 

A latent assumption in most previous works is that ev- 
ery node in a transport network plays the role of a host, 
that is to say, every node has the ability creating a certain 
kind of substance, energy or information. However, the 
real world is far from this assumption. For example, in 
an electric network [HI, [13], there are two kinds of nodes, 
power stations and transformer substations. The power is 
generated in the former nodes, flowing to the latter ones, 
and then imported to customers through them. There- 
fore, power stations behave as a kind of suppliers, while 
the transformer substations are customers holding de- 
mands. In some Internet serving systems, such as music 
libraries (e.g., audioscrobbler.com, see Ref. [i3|), movie- 
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sharing systems (e.g., Netflix.com, see Ref. ^1^]) and 
on-line viewing site (e.g., YouTube.com, see Ref. [isj). 
all the resources are located in a few servers, while other 
connected nodes, usually personal computers, only regale 
themselves with those services. Those examples give rise 
to a general concept of supply- demand network, whose 
nodes are classifled into two categories: the supply nodes 
provide some kinds of services, while the demand nodes 
play the role of customers. Analysis of supply-demand 
networks has found its applications in various real sys- 
tems, ran ging from the power grid [l6l.[l7| to supply chain 
networks |l8l . [Tgj . 

In this paper, we propose a general model for the trans- 
port on a supply-demand network, whose capacity is very 
sensitive to the locations of supply nodes. By applying a 
simulated annealing algorithm, we obtained the near op- 
timal locations of supply nodes subject to the maximal 
network transport capacity. The proposed algorithm per- 
forms obviously better than the random selection, degree 
targeted, betweenness targeted, and greedy methods. 

II. MODEL 

Considering a network consisted of N nodes, which are 
classified into two categories: One is called the supplier 
that provides a certain kind of service, the other is called 
the customer who requires this service. Here, the ser- 
vice is an abstract concept and can stand for substance, 
energy, information, etc. For simplicity, we use the lan- 
guage of the Internet, that is to say, every customer need 
some information packets (resource) , and only the suppli- 
ers can generate those information packets. We assume 
the demands are uniformly distributed, namely each cus- 
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FIG. 1: Illustration of the distribution of edge loads in a 
supply-demand network. The gray solid and hollow circles 
denote supply nodes and demand nodes, respectively. In each 
of panels (a), (b) and (c), the circle marked by a star is the 
target demand node, and the resulting loads are labeled be- 
sides corresponding edges. Integrating (a), (b) and (c), the 
distribution of edge loads can be obtained, as shown in the 
panel (d). Here, Lmax ~ 4/3. 

tomer needs a unit number of packets (one can simply 
say one packet). For a given customer, we suppose this 
packet is always sent by one of the nearest suppliers. 
However, in general case, there are several nearest suppli- 
ers and for each there are several shortest paths. In the 
real implementation, one of those shortest paths should 
be randomly picked, and the packet will follow this path 
from the supplier to the customer. In the numerical cal- 
culation, to reduce the fluctuation, if there are in parallel 
k shortest paths from a customer to the suppliers (gener- 
ally, those paths aim to more than one nearest suppliers), 
we assume the packet is divided into k pieces, each goes 
through one shortest path and contributes l/Zc to the 
traffic load (see an illustration in Fig. 1). 

If the bandwidth (i.e., traffic capacity) of each edge is 
identical, the maximal edge load, Lmax, is the key factor 
determining the traffic condition. Actually, the traffic 
congestion will occur when Lmax exceeds the bandwidth. 
Therefore, given a limited bandwidth, the smaller Lmax 
corresponds to higher transportation capacity. Analo- 
gously, in the previous studies [S, 7], the maximal node 
load is usually used to quantify the system's performance: 
the smaller the maximal node load, the higher the trans- 
port capacity. In this paper, we use edge load instead 
of node load because in the real systems, such as the In- 
ternet and the highway, the congestion usually happens 
along the edges, not at the nodes [20j . 

Given a network structure and the number of suppliers, 
we aim at finding out the optimal configuration of suppli- 
ers (i.e., the locations of suppliers) making L^ax as small 
as possible. This is an optimization problem (indeed, an 
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FIG. 2: The objective function, Lmax, vs. system time, t, in the 
optimizing process of the SA algorithm. This figure illustrates 
a typical result on a BA network of size N = 1000, average 
degree (fc) — 6. The number of suppliers is set as M = 10. 

NP hard problem) with Lmax being the objective function, 
and the algorithm presented in this paper (see below) can 
be directly extended to the case with maximal node load 
being the objective function. In addition, since many 
real transportation networks have heterogeneous degree 
distribution (see the examples shown in Refs. [2]]. |22|). 
we use scale-free networks to mimic their topologies. 



III. ALGORITHM 

In a supply-demand network of N nodes and M sup- 
f N\ . 

pliers, there are in total ( I different configurations 

for suppliers' locations. Finding the optimal solution 
by evaluating all the possible configurations is infeasible 
when TV 3> M ^ 1. The optimization of a system with 
many degrees of freedom with respect to a certain objec- 
tive function is a frequently encountered task in physics 
and beyond. One special class of algorithms used for find- 
ing the high-quality solutions to those NP-hard optimiza- 
tion problems is the so-called nature inspired algorithms, 
including simulated annealing (SA) [23l |2^. genetic al- 
gorithms (GA) [25„ ,26.], genetic programming (GP) [27| . 
extremal optimization (EO) [2R |29||. and so on. Here we 
adopt the SA algorithm, whose procedure is as follows. 

(i) Randomly choose an initial configuration, denoted 
by S^. Calculate its maximal edge load, i^axi ^^'^ set 
the best solution as: S^^^"" = S° and L^^f = L^^^. Set 
the system time as i = 1. 

(ii) Randomly pick one supplier from the configuration 
S*~^, and change its location randomly, denote this new 
configuration as S''. Calculate its maximal edge load, 

Tt 
max" 

(ui) If < Llll\ then set 5^-^ = and LlH' = 
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TABLE I; Comparison of the maximal edge load obtained by 
DTA, BTA, GM and SA. The underlying networks are BA 
networks with TV — 1000 and (fc) — 6, and all the data are 
obtained by averaging over 100 network configurations. 



M / Algorithm DTA BTA GM SA 
5 14.73 14.98 13.32 12.37 

10 8.25 8.92 7.17 6.31 



-^max- If -^max ^ ^iax ' accept the Current configura- 
tion, that is, set t ^ i + 1 and repeat (ii). Otherwise, if 
> Ll^^, the current configuration is accepted with 
probabiUty e~^/^, where T is a temperature-like param- 
eter and A — L^^^^ — L^^^. When a configuration is re- 
jected, the algorithm directly goes back to (ii) and keeps 
the system time t unchanged. 

To obtain the high-quality solution, one shall repeat 
the step (ii) as long as desired. In this paper, we termi- 
nate the algorithm if the variance of i^ax ™ the latest 
10^ time steps is smaller than a threshold 10^^. Note 
that, one time step corresponds to one implementation 
of step (ii), which is different from the system time t. 
The parameter T is crucial for the algorithmic efficiency. 
According to the MatropoHs's Guidance [13], in the ini- 
tial stage, the accepting probability of a new configura- 
tion should be close to 1. Therefore, we first choose a 
relatively low temperature Tq, and numerically calculate 
the corresponding accepting probability, resulted from a 
random change of one supplier's location in a completely 
random configuration. The temperature is doubled un- 
til the accepting probability reaches a threshold quan- 
tile 0.50. During the searching process, the temperature 
should slowly decrease [l^l, here we adopt the simplest 
method, that is, we set T ^ aT after every Q time steps, 
where the parameter a is 0.90 and the period is set as 
Q = O.liVM. 

For comparison, we also implement some other algo- 
rithms. A brief introduction is as follows. Random Al- 
location (RA) - The locations of suppliers are selected 
completely randomly. Degree Target Algorithm (DTA) - 
The suppliers are set as the M nodes with highest de- 
grees. Betweenness Target Algorithm (BTA) - The sup- 
pliers are set as the M nodes with highest betweennesses 
(see Refs. [3l|, HI] for the definition and calculation of 
node betweenness). Greedy Method (GM) - First, we 
consider the case with only one supplier, and find out 
the optimal location of this supplier that minimizes the 
corresponding L^ax- Then, we add one supplier and find 
out its optimal location under the condition that the lo- 
cation of the firstly added supplier is fixed. Repeating 
this operation, that is, at the fcth step, we add the kth 
supplier and find out its optimal location subject to mini- 
mal Lmax under the condition that the locations of former 
k — 1 suppliers are fixed. This algorithm is terminated 
when M suppliers are added already. 
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FIG. 3: (Color online) Algorithmic performance for BA net- 
works. The main plot shows a comparison among DTA, BTA, 
GM and SA, while the inset reports a comparison between 
RA and SA. The number of suppliers, M, varies from 1 to 
10, while the network size = 1000 and the average degree 
(fc) = 6 are fixed. All the data points are obtained by aver- 
aging over 100 network configurations. 
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FIG. 4: (Color online) Scatter plot of betweenness vs. degree 
in a BA network with N — 1000 and (k) = 6. Each small 
black fork represents a node. These 10 red circles denote the 
selected suppliers by SA. The smallest degree of suppliers is 
9, and the second smallest one is 12. 



IV. RESULTS 

In this paper, all the numerical simulations are imple- 
mented based on the Barabasi- Albert (BA) model [33| . 
which is one of the minimal models reproducing the het- 
erogenous structure of real- world networks. Figure 2 
reports a typical optimizing process, during which the 
objective function, Lmax, fluctuates strongly in the early 
stage and approaches to a relatively stable value lately. 
The proposed SA can reduce the objective function, Lmax, 
by more than 10 times compared with its initial value 
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corresponding to a random selection of suppliers. We 
implement SA in larger BA networks {N = 1000) for dif- 
ferent M from 1 to 10, and take the average over 100 
independent network configurations. As shown in the in- 
set of Fig. 3, SA performs much better than RA. We 
also compare SA with some mentioned algorithms, DTA, 
BTA and GM, and the results have demonstrated that 
SA performs best. We report two examples, M = 5 and 
M = 10, in Table 1. The improvement is in general 
about 10%. Note that, although SA performs the best, 
it spends the longest running time. Actually, the time 
complexity obeys the inequality 0{SA) > 0{GM) > 
0{BTA) > 0{DTA). Since GM performs not so bad, 
it is a strong candidate especially for huge-size networks, 
and GM might be a considerable tradeoff of time com- 
plexity and accuracy of solution. 

Note that, although BA model has successfully cap- 
tured the degree heterogeneity of real networks, it lacks 
some other important structural properties, such as 
the community structure [s^l and rich-club phenomenon 
[sH - DTA might perform worse if the network has 
strongly community structure or presents the rich-club 
phenomenon. The reason is a good algorithm should pre- 
fer to allocate suppliers to different communities rather 
than putting them together in a community contain- 
ing many very-large-degree nodes, and if the very-large- 
degree nodes are closely connected to form a rich club, 
selecting them as a whole is of low efficiency since the 
increasing suppliers cannot substantially reduce the av- 
erage distance from customers to suppliers. As a start 
point, we here only discuss simulation results on BA net- 
works, and leave the investigations of algorithmic perfor- 
mance on more complicated topologies as an open issue. 

The DTA and BTA have almost the same performance 
and give out very similar selections of suppliers, for in 
BA networks betweenness and degree are very strongly 
correlated [s^ [s^l • To provide insights of the solution by 
SA, in Fig. 4, we give a scatter plot of betweenness versus 
degree, and mark by red those selected suppliers. Though 
SA also prefers large-degree (large-betweenness) nodes, 
the selected suppliers are remarkably different from those 
by DTA or BTA, actually, moderate-degree (moderate- 
betweenness) nodes also have chance to be selected by 
SA. In most cases, only the top-40% large-degree nodes 
have the chance to be selected, therefore we can restrict 
the candidates of suppliers in those 40% nodes. We have 
check this restriction in BA networks with N = 1000 and 
{k) = 6, which gives out equivalently good solution while 
requires about 10 times shorter CPU time. 



V. CONCLUSION AND DISCUSSION 

In this paper, we proposed a generic model of transport 
in supply-demand network, which is consisted of suppli- 
ers (supply nodes) and customers (demand nodes). Ac- 
cordingly, a measure of edge load is given, under the as- 
sumption that every customer only requires service from 



the nearest supplier. In such a network with heteroge- 
nous degree distribution, its transport capacity is very 
sensitive to the locations of supply nodes. We there- 
fore design a simulated annealing algorithm to find out 
the near optimal configuration of supply nodes, which 
remarkably enhances the transport capacity, and outper- 
forms the degree target algorithm, the betweenness tar- 
get algorithm, and the greedy method. This work pro- 
vides a start point for systematically analyzing and opti- 
mizing transport dynamics on supply-demand networks. 
Even though the model and algorithm arc simple, we get 
some non-trivial result, that is, simply picking up those 
nodes of highest degrees is not the optimal method, ac- 
tually, some moderate-degree nodes also have chance to 
be selected as suppliers. 

In our model, every customer requires the same 
amount of resource, which is not in accordance with the 
elephants and mice phenomenon ^SS'l found in the real 
Internet, where a small fraction of flows contribute to 
most of the traffic. Corresponding to the current model, 
a flow stands for the resource transported from a supplier 
to a customer, and thus each flow has the same size al- 
though the one passing longer paths contributes more to 
the total load. In addition, the proposed algorithm does 
not fully take into account and make use of the topo- 
logical features. We have already mentioned in the last 
section that the mesoscopic structure, such as communi- 
ties and the rich club, may highly influence the solutions. 
Those structural information should be extracted prior to 
the optimizing algorithm, and be embedded in the algo- 
rithmic procedure in some way to improve the efficiency 
and/or the resulting network capacity. All those blem- 
ishes listed above can be treated as some open problems 
worth of a future exploration. 

To the end, we emphasize that many real systems can 
be better described by the current supple-demand net- 
work model, instead of the oversimple assumption Q 
that every node simultaneously plays the roles of sup- 
plier and customer. We have already mentioned some 
examples, such as power grid [l^, [l3| and supply chain 
networks [H, [l^ , another typical example is the software 
supporting systems in the Internet, where a system usu- 
ally has set up several servers in different locations, and 
users from everywhere can ask for downloading of some 
softwares. The locations of those servers play the cru- 
cial role in determining the efficiency and capacity the 
software supporting system. 

This study also provides some complementary infor- 
mation for relevant phenomena in disparate systems. For 
example, social scientists have studies how to design who 
should be integrators in a given social communication 
networks to better solve problems, and they have found 
that people having extensive relations (i.e., of very large 
degrees) may not be the suitable information integrators, 
instead, the highest efficient structure makes the distance 
of all nodes from the obvious integrator the shortest [s^ , 
which is, to some extent, in accordance with what we ob- 
served in this work. In addition, empirical studies show 
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that the pubhc service facilities are not just located in the 
place of the most dense population, but somehow more 
uniformly distributed to make the total travel distance 
between people and facilities shorter . As a final remark, 
we noted that a very recent work has considered of the 
network-based transport with multiple sources and sinks 
[4l| , which shows different yet relevant motivation to the 
current work. 
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